Phase behavior and far-from-equihbrium gelation of charged attractive colloids 
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In this Rapid Communication we demonstrate the apphcabihty of an augmented Gibbs ensemble 
Monte Carlo approach for the phase behavior determination of model colloidal systems with short- 
ranged depletion attraction and long-ranged repulsion. This technique allows for a quantitative 
determination of the phase boundaries and ground states in such systems. We demonstrate that 
gelation may occur in systems of this type as the result of arrested microphase separation, even 
when the equilibrium state of the system is characterized by compact microphase structures. 
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Understanding the routes to colloidal gel formation 
is a problem in the forefront of soft-condensed matter 
physics. It has been demonstrated via theory [1], experi- 
ment [2|,|3[, and simulation y,|5[ that a weak, porous solid 
may be formed when a suspension of colloidal particles 
with short-ranged attractions is quenched below its crit- 
ical point. Spinodal decomposition, expected to lead to 
complete phase separation into colloid-rich and colloid- 
poor regions, may arrest or become anomalously slow 
when bonding between colloidal particles is sufficiently 
strong. The resulting structure may support weak shear 
stresses, and is called a colloidal gel. While equilibrium 
routes to physical colloidal gelation do exist [y, [Zl, |8|, 
the nonequilibrium route described above is perhaps the 
most ubiquitous. 

When excess charge resides on the colloidal particles, 
the situation is more complex [9[. In addition to short- 
ranged depletion interactions, the long-ranged repulsive 
portion of the potential may lead to the formation of var- 
ious microphase structures, such as clusters, cylinders, 
sheets, and spirals [lO|, UjJ, ll2|, [isl, [l^]- As a result, the 
nature of physical gelation in these systems is potentially 
more subtle than in systems with attractive interactions 
alone. Indeed, several novel routes to gelation in these 
systems have been discussed. One proposal is that com- 
pact, thermodynamically stable clusters may form the 
building blocks of a glassy state [l5[ . At high enough vol- 
ume fraction clusters may aggregate or become caged, as 
in the case for individual colloidal particles near the col- 
loidal glass transition [l6[. At low volume fractions, such 
a glassy state could be stabilized by the effective long- 
ranged repulsion between individual clusters [l5[. More 
recently, it has been proposed that disordered or partially 
ordered anisotropic domains may be responsible for phys- 



ical gelation in charged colloidal systems |llUl4lJl7[. 

A major difficulty in testing the validity of any pro- 
posal resides in the accurate calculation of the phase di- 
agram in systems with long-ranged repulsive and short- 
ranged attractive interactions. When the range of the 
attractive interaction is comparable to the size of the 
particles, the phase behavior may be characterized eas- 
ily with standard simulation techniques [12, |l3|, [jj. On 
the other hand, when the attractive portion of the po- 
tential varies over distances that are a fraction of the 
particle size, sampling issues become severe, and direct 
molecular dynamics (MD) and Monte Carlo (MC) tech- 
niques are not feasible. Unfortunately it is in the short- 
ranged attractive limit that the gel phase is usually stabi- 
lized [4, _5]. We have recently demonstrated [18, 19] that 
a judicious implementation of the Gibbs ensemble Monte 
Carlo (GEMC) [20] technique allows for phase behavior 
computation in systems with very short-ranged attrac- 
tive interactions. In this Rapid Communication we adapt 
the GEMC technique for the study of systems character- 
ized by short-ranged attraction and long-ranged repul- 
sion. The successful application of this approach provides 
crucial information on possible dynamical routes to gela- 
tion in systems with competing interactions, which are 
examined via MD simulation. Indeed, our conclusions 
and interpretations of the gelation mechanism differ from 
those previously drawn from studies of analogous sys- 
tems. 

In this work we study potentials of the form^ ^5[ 



U{r) = 4e 



iT-{^)' 



r/i 



r/i ' 



(1) 



where large n values correspond to short-ranged attrac- 
tions. We study three-dimensional systems (n=50) with 
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^ Temperature T is in units of well depth e, distance in particle 
diameter cr, and time t in ^/e/ma'^ with particle mass m. 
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FIG. 1: (Color online) Phase diagrams for short-ranged at- 
traction in (a) three (n=50) and (b) two (n=100) dimensions 
with increasing repulsion set by A and ^. The dashed line 
is the crystal close-packed density. Critical points (open cir- 
cles) for purely attractive systems (A=0) are obtained as de- 
scribed in the text. Small symbols indicate gas-solid and gas- 
Uquid coexistence and large symbols the microphase bound- 
ary for ^=2 and A=0.2. Uncertainty is smaller than the sym- 
bols, unless otherwise indicated. The phases are labeled M 
(microphase), H (homogeneous fluid), C (gas-crystal coexis- 
tence), and S (solid). Lines and shaded areas are guides for 
the eye. Inset: F^/c, t) = (exp {k • [r^(t) - ri{0)]}) from MD 
of the microphase forming system at the \ow-k peak of the 
structure factor ka^O.9 for p=0.6 and T=0.22, 0.17, 0.12, 
and 0.08, from left to right. For the higher two temperatures, 
trajectories are run for many times the relaxation time. For 
the lower two temperatures the system shows aging. 



three sets of repulsive Yukawa parameters {A=0 [l8[; 
A=1.275, C=0.5; A=0.203, ^=1) and two-dimensional 
systems (n=100) with parameters (A=0; 74=0.68, ^=1; 
74=0.2, ^=2). For the last parameter set, the system is 
beyond the Lifshiftz point and thus microphase separates 
at low temperatures. A preliminary dynamical study of 
a related three-dimensional system has been performed 
by Sciortino et al. ^]. For visualization purposes, we 
focus on a two-dimensional system. As discussed below 
we expect all qualitative conclusions to be unmodified by 
this choice. Selected results will also be provided for the 



three-dimensional case. 

For the vapor-solid equilibrium, we use the GEMC 
methodology described in Ref.^ [18]. Initial boxes con- 
tain 256 particles in the gas phase and 500 particles in 
the solid slab. A minimum of 10^ MC cycles are per- 
formed for equilibration and production. To characterize 
the phase behavior in the regime where microphase sepa- 
ration occurs, two identical boxes with a minimum of 256 
particles on a lattice are used for initial configurations. 
At these lower temperatures, the GEMC method serves 
as an efficient phase space sampling algorithm akin in 
spirit to the grand canonical MC and parallel-tempering 
MC used in previous studies of systems with compet- 
ing attractive and repulsive interactions [13]. It should 
be noted, however, that the approach implemented here 
is more efficient at detecting complex phase boundaries. 
In particular, upon crossing from a gas-solid coexistence 
region to a microphase separated region, the parallel- 
tempering MC sampling efficiency is dramatically re- 
duced. 

Before discussing the impact on phase equilibrium of 
the repulsive portion of the potential, we make some 
brief comments on the purely attractive case in two di- 
mensions. The phase diagrams of short-ranged attrac- 
tive two-dimensional systems have not been investigated 
as thoroughly as their three-dimensional counterparts. 
The existence of high-density crystal-crystal coexistence 
in both two and three dimensions supports the notion 
that systems characterized by identical interactions in 
two and three dimensions are qualitatively similar [22[. 
While crystal-crystal coexistence does not occur in the 
temperature and density range of interest here, this phe- 
nomena is related to the metastability of the gas-liquid 
critical point, which does occur at densities and tem- 
peratures of concern in this work. With no repulsion 
(^4=0), we thus expect the liquid-gas coexistence to be- 
come metastable for n > 12, so for n=50-100 the gas- 
liquid binodal is metastable and buried below the gas- 
solid coexistence line [23]. In three dimensions it is known 
that nucleation dynamics are dramatically influenced by 
this buried gas-liquid binodal [2J, |25|, |26[. Similar ef- 
fects also appear in two dimensions, although metastable 
states are harder to observe, due to much faster nucle- 
ation. In the two-dimensional attractive case (n=100, 
^4=0), our visual inspection of the relaxation dynam- 
ics accompanying apparent spinodal decomposition gives 



^ Local and non-local moves are equally frequent at low densities, 
high temperatures, and for gas-liquid equilibria. Otherwise, non- 
local moves comprise up to 90% of the total. Local moves are 
particle displacements (94%), cluster displacements (0.5%), rota- 
tions (0.5%), and cleaving |3l| (up to 1%), as well as symmetric 
(3%) and asymmetric (1 to 2%) volume exchanges. Non-local 
moves are particle exchanges (80%) and one or two-box (10% 
each) aggregation volume bias [3^. 
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FIG. 2: (Color online) Configurations (^=2, A=0.2) equilibrated with GEMC simulations at (a) p=0.3, (b) p=0.5, and (c) 
y9=0.6 as well as (d) p=0.6 MD evolution after quench to T=0.12 (first two panels) and GEMC thermodynamic outcome (third 
panel) . 



Tc=0.18±0.01. This can be favorably checked by combin- 
ing the arguments developed by Noro and Frenkel regard- 
ing corresponding states for short-ranged potentials [2l| 
with the location of the critical point coordinates ex- 
tracted from the results of Seaton and Glandt for the 
Baxter limit of adhesive disks [27], which further give 
an estimate of the critical density pc=0.6 ±0.1. For the 
same potential form (n=100, A=0), but in three dimen- 
sions, Tc=0.235 ± 0.005 [15] and pc^O.b 0, [1^, while 
the three-dimensional system studied here (n=50, A=0) 
has Tc=0.288± 0.005 Q- The metastable critical points 
are indicated in Fig. [H 

We now turn to a discussion of the impact of repul- 
sion on phase behavior. In three-dimensional systems 
with attractive ranges similar to that of of the n=6 
(Lennard- Jones) potential, the addition of small ampli- 
tude long-ranged repulsion is known to slightly depress 
the gas-liquid spinodal and binodal lines [29]. Repul- 
sion of greater amplitude or longer range transforms the 
depressed gas-liquid coexistence line into a microphase 
boundary where low- and high-density regions locally 
coexist, as in Fig. W[h). Besides, repulsion destabilizes 
the crystal phase and thereby extends the stability of 
the homogeneous fluid phase to lower temperatures, as 
seen in Fig. [H Where microphase separation occurs, this 
depresses the low-density microphase boundary. In this 
case the high-density phase is crystalline, since the mi- 
crophase boundary is below the gas-solid coexistence line. 
GEMC temperature scans allow us to visually pinpoint 
the boundary of the microphase region. General argu- 
ments indicate that fluctuations alter the nature of the 
critical point, rendering the transition weakly first or- 
der [3^. 

At this stage it is useful to directly observe representa- 
tive equilibrium configurations obtained during GEMC 
runs, as in Figs. [2l[a)-(c). Configurations below the mi- 
crophase transition show slightly irregular domains, due 








FIG. 3: (Color online) For p=0.5 (^=2, A=0.2) (a) a nearly 
arrested gel configuration obtained by MD at t=5 x 10^ after 
a quench to T=0.08 and (b) an equilibrium configuration at 
T=0.08 found by GEMC simulation. 



to finite-temperature entropic fluctuations. For this sys- 
tem at the densities and temperatures shown, compact 
cluster states are more stable than extended (lamellar) 
states, which were not observed. The second and third 
panels of Fig. ^ illustrating configurations just above 
and below the microphase line are particularly illuminat- 
ing. The dynamical behavior upon a rapid quench below 
the microphase separation line is shown in the first two 
panels of Fig. EJd).^ The evolution of the system from 
the second to the third panel of Fig. Efd) is extremely 
slow via conventional MD. The dramatic increase in time 
scales upon cooling the system below the microphase sep- 
aration can be seen from the self-intermediate scattering 



^ MD integration with step At=0.001 was performed for systems 
sizes A^ > 256. Cooling is done from a fluid configuration to 
the set temperature by resampling velocities every tcool=^^- 
The type of dynamics does not influence the qualitative results 
strongly [5], so other protocols are not studied. Slower cooling 
rates allow the system to reach configurations closer to equilib- 
rium [l^. The potential is truncated at 7-0=5(3.3) for ^=2(1), 
which alters the microphase patterning J13l| , but here the results 
do not depend on the form of the ground state, since GEMC and 
MD runs use equal Tc- 



function Fs{k^t) in the inset of Fig. [T](b). Fluctuations 
on the order of the characteristic domain length scale de- 
cay so slowly that the structure is effectively dynamically 
arrested. 

The dynamical behavior shown in Fig. Efd) is clearly 
relevant for the mechanism of gelation of the system. In 
particular, consider slightly deeper quench at a lower den- 
sity than that shown in Fig. WiA)- In Fig. [Sja) a quench 
of the system is made at p=0.5 to T=0.08. After rapid 
initial transient dynamics, where a tenuous percolating 
structure is formed, the dynamics becomes anomalously 
slow. This system may be characterized as a gel. Due to 
the depth of the quench, we did not observe substantial 
structural coarsening [as is observed in Fig. El^d)] over 
the duration of the simulations, although the energy of 
the system is slowly evolving in time. The similar evo- 
lution is expected to be much slower in the case of a 
polydisperse, three-dimensional system. For comparison, 
the equilibrium state of the system (as found via GEMC 
simulation) is also shown in Fig. [3l Note that the struc- 
ture of the gel does not resemble the equilibrium state 
in any regard, despite the fact that the energy per par- 
ticle in the gel is quite close to that of particles in their 
cluster ground state. In addition, we find via the GEMC 
method no evidence of extended microphase structures 
at the density and temperature of Fig. [31 Regardless, 
the short-time critical-like fluctuations that exist during 
the initial stages of microphase separation are sufficient 
to generate configurations that span space and evolve 
anomalously slowly, due to strong, short-ranged bonding 
between particles. 

The results presented above are clearly of importance 
for understanding the general routes to gelation in col- 
loidal systems with competing short-ranged attraction 
and long-ranged repulsion. We find, rather robustly, that 
gels may form via arrested microphase separation. This 
is analogous to the purely attractive case, where gels are 
formed by arrest of global phase separation. This picture 
is in contrast to that of Refs. [l5|, ll6|, where gelation is 
a consequence of near-equilibrium vitrification of clusters 
that are stabilized by the repulsive interactions. We have 
explicitly demonstrated here that, as in the purely at- 
tractive case, the gel may be stabilized by the attraction 
between particles. More interestingly, the gel structure 
bears no similarity to the clusters that exist in thermo- 
dynamic equilibrium. We find that in the system stud- 
ied here, the only role of the repulsion is to select the 
symmetry of the ground state(s). Our results are fully 
consistent with, but more general than, the suggestion 
that disordered configurations of anisotromc microphase 
structures are connected to gelation Uj, |lj, ll7| . Indeed, 
we have directly demonstrated that extended microphase 
structures are not involved during the formation of the 
gel in Fig. [3l Furthermore, preliminary studies in three 
dimensions for the system studied in Ref. [15] show that 
for rapid, deep quenches gels may be formed in the same 



manner. This leads us to believe that both the routes 
presented here and in Refs. [15', 'iGl] are possible, and will 
likely depend intimately on the details of the quench pro- 
tocol. In general, lamellar and tubular structures (as well 
as other microphase textures) may comprise the building 
blocks of the disordered weak solids found in colloidal sys- 
tems with depletion attraction and long-ranged charge 
repulsion, but our results demonstrate that they need 
not. The disordered lamellar and tubular examples are 
simply some of the many ways in which microphase sep- 
aration may be arrested in a manner directly analogous 
to the arrest of global phase separation. The collection of 
all such ways perhaps forms the simplest generic route to 
the gel state at intermediate volume fractions in attrac- 
tive systems with charge frustration. The search for more 
exotic routes which may be operative at very low volume 
fractions and where charge repulsion plays a more active 
role may be greatly facilitated by the use of the GEMC 
approach as outlined here. 
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